Origen de los datos

Los datos son provenientes del gobierno de Mexico para el estudio de aguas subterraneas, en el rango de tiempo del 2012-2024. Muestra los distintos indicadores que encuentran en el agua

Los Indicadores de calidad del agua son una herramienta cuantitativa que utiliza la Gerencia de Calidad de Agua de la Conagua para determinar cómo se encuentra la calidad del agua en diversos sitios de los cuerpos de agua nacionales, clasificados como lóticos, lénticos, costeros o subterráneos. Se determinan a partir de los resultados del monitoreo nacional que lleva a cabo la Red Nacional de Medición de la Calidad del Agua (RENAMECA), operada por la Gerencia de Calidad del Agua, de la Subdirección General Técnica.

Los Indicadores subterráneos son 14 parámetros: Fluoruros (FLUORUROS), Coliformes fecales (COLI_FEC), Nitrógeno de Nitratos (N_NO3), Arsénico Total (AS_TOT), Cadmio Total (CD_TOT), Cromo Total (CR_TOT), Mercurio Total (HG_TOT), Plomo Total (PB_TOT), Alcalinidad (ALC), Conductividad eléctrica (CONDUCT), Dureza Total (DUR), Sólidos Disueltos Totales (SDT), Hierro Total (FE_TOT) y Manganeso Total (MN_TOT).

Asimismo, se utiliza un Semáforo de la calidad del agua, que es una clasificación de la calidad en un determinado sitio de monitoreo, basada en intervalos de concentración de cada parámetro indicador, por lo que se establece una escala de cumplimiento que asigna un color: verde que representa la mejor calidad, amarillo calidad intermedia y rojo mala calidad.

El semáforo para agua subterránea considera lo siguiente: si uno o más de los siguientes parámetros: FLUORUROS, COLI_FEC, N_NO3, AS_TOT, CD_TOT, CR_TOT, HG_TOT y PB_TOT, incumplen, entonces el semáforo será color rojo; si los parámetros anteriores cumplen, pero uno o más de los siguientes parámetros: ALC, CONDUCT, DUR, SDT, FE_TOT y MN_TOT, incumplen, el semáforo será amarillo; y, cuando se da el cumplimiento de todos los indicadores, el semáforo será verde.

Analisis de Datos

1 - Cargamos los Datos

Debido a que el archivo con los datos se encuentra en formato xlsx utilizamos la biblioteca que nos permite leerla

library(readxl) #Paquete para leer excel
library(tidyverse) #Libreria principal para el analisis de datos
library(ggplot2)
library(scales)
library(plotly)
library(ggcorrplot)
library(randomForest)
library(caret)
library(class)      # Para la función knn
library(caret)      # Para partición y evaluación
data_agua <- read_excel("Calidad_del_Agua_Subterr_nea_p_2012-2024_15082025.xlsx")

2 - Analisis de Datos Crudos

Realizando una inspeccion a los datos que trae el archivo vemos que trae 41 columnas donde: - Las primeras 6 columnas traen informacion sobre el lugar donde se registraron esos datos - Las siguientes columnas son distintos indicadores quimicos/biologicos para determinar la calidad del agua en esos lugares - Al lado de cada indicador hay una etiqueta que clasifica ese indicador particular para determinar su influencia en la calidad del agua - La ultimas dos columnas son un semaforo que clasifica en VERDE, AMARILLO, ROJO las muestras de agua dependiendo del total de los indicadores

head(data_agua)
## # A tibble: 6 × 41
##   `CLAVE SITIO`   SITIO `ORGANISMO _DE_CUENCA` ESTADO MUNICIPIO ACUIFERO SUBTIPO
##   <chr>           <chr> <chr>                  <chr>  <chr>     <chr>    <chr>  
## 1 BROTE CARMINA 3 CARM… RÍO BRAVO              COAHU… ACUÑA     PRESA L… BROTE  
## 2 CAZEPA-1        POZO… CUENCAS CENTRALES DEL… DURAN… LERDO     VILLA J… POZO   
## 3 CAZEPA-10       POZO… CUENCAS CENTRALES DEL… DURAN… LERDO     PRINCIP… POZO   
## 4 CAZEPA-11       POZO… CUENCAS CENTRALES DEL… DURAN… LERDO     PRINCIP… POZO   
## 5 CAZEPA-2        POZO… CUENCAS CENTRALES DEL… DURAN… LERDO     VILLA J… POZO   
## 6 CAZEPA-25       POZO… CUENCAS CENTRALES DEL… DURAN… LERDO     VILLA J… POZO   
## # ℹ 34 more variables: LONGITUD <dbl>, LATITUD <dbl>, PERIODO <chr>,
## #   `ALC_mg/L` <dbl>, CALIDAD_ALC <chr>, `CONDUCT_mS/cm` <dbl>,
## #   CALIDAD_CONDUC <chr>, `SDT_mg/L` <dbl>, CALIDAD_SDT_ra <chr>,
## #   CALIDAD_SDT_salin <chr>, `FLUORUROS_mg/L` <dbl>, CALIDAD_FLUO <chr>,
## #   `DUR_mg/L` <dbl>, CALIDAD_DUR <chr>, `COLI_FEC_NMP/100_mL` <chr>,
## #   CALIDAD_COLI_FEC <chr>, `N_NO3_mg/L` <dbl>, CALIDAD_N_NO3 <chr>,
## #   `AS_TOT_mg/L` <dbl>, CALIDAD_AS <chr>, `CD_TOT_mg/L` <dbl>, …

3 - Procesamiento de Datos

Vamos a eliminar algunas columnas que no aportan al estudio que se esta realizando. Para no complejizar este analisis eliminamos las etiquetas de cada uno de los indicadores, dejando solo el semaforo final

Tambien eliminamos algunas de las columnas con informacion del lugar donde se tomo la muestra, dejando solo las mas generales para poder hacer un analisis simplificado

Luego renombramos las columnas para mejorar la visualizacion en los diferentes graficos

datos_limpios1 <- data_agua %>%
  select(-c("CLAVE SITIO","SITIO","MUNICIPIO","LONGITUD","LATITUD","PERIODO","CONTAMINANTES","CALIDAD_ALC","CALIDAD_CONDUC","CALIDAD_SDT_ra","CALIDAD_SDT_salin","CALIDAD_FLUO","CALIDAD_DUR","CALIDAD_COLI_FEC","CALIDAD_N_NO3","CALIDAD_AS","CALIDAD_CD","CALIDAD_CR","CALIDAD_HG","CALIDAD_PB","CALIDAD_MN","CALIDAD_FE"))
datos_limpios2 <- datos_limpios1 %>%
  rename("ALC"="ALC_mg/L","CONDUCT"="CONDUCT_mS/cm","SDT"="SDT_mg/L","FLUORUROS"="FLUORUROS_mg/L","DUR"="DUR_mg/L","COLI_FEC"="COLI_FEC_NMP/100_mL","N_NO3"="N_NO3_mg/L","AS"="AS_TOT_mg/L","CD"="CD_TOT_mg/L","CR"="CR_TOT_mg/L","HG"="HG_TOT_mg/L","PB"="PB_TOT_mg/L","MN"="MN_TOT_mg/L","FE"="FE_TOT_mg/L","SEMAFORO"="SEMÁFORO")
head(datos_limpios2)
## # A tibble: 6 × 19
##   `ORGANISMO _DE_CUENCA`   ESTADO ACUIFERO SUBTIPO   ALC CONDUCT   SDT FLUORUROS
##   <chr>                    <chr>  <chr>    <chr>   <dbl>   <dbl> <dbl>     <dbl>
## 1 RÍO BRAVO                COAHU… PRESA L… BROTE    140.     884  562.     0.82 
## 2 CUENCAS CENTRALES DEL N… DURAN… VILLA J… POZO     149.     329  185.     0.52 
## 3 CUENCAS CENTRALES DEL N… DURAN… PRINCIP… POZO     205.     615  378      0.43 
## 4 CUENCAS CENTRALES DEL N… DURAN… PRINCIP… POZO     229.     602  322      0.5  
## 5 CUENCAS CENTRALES DEL N… DURAN… VILLA J… POZO     168.     379  240.     0.566
## 6 CUENCAS CENTRALES DEL N… DURAN… VILLA J… POZO     170.     354  234.     0.48 
## # ℹ 11 more variables: DUR <dbl>, COLI_FEC <chr>, N_NO3 <dbl>, AS <dbl>,
## #   CD <dbl>, CR <dbl>, HG <dbl>, PB <dbl>, MN <dbl>, FE <dbl>, SEMAFORO <chr>
# Ver un resumen estadístico de todas las columnas
summary(datos_limpios2)
##  ORGANISMO _DE_CUENCA    ESTADO            ACUIFERO           SUBTIPO         
##  Length:2728          Length:2728        Length:2728        Length:2728       
##  Class :character     Class :character   Class :character   Class :character  
##  Mode  :character     Mode  :character   Mode  :character   Mode  :character  
##                                                                               
##                                                                               
##                                                                               
##                                                                               
##       ALC             CONDUCT              SDT            FLUORUROS      
##  Min.   :  10.79   Min.   :   11.13   Min.   :    6.0   Min.   : 0.0300  
##  1st Qu.: 155.90   1st Qu.:  442.00   1st Qu.:  310.0   1st Qu.: 0.3129  
##  Median : 211.81   Median :  721.00   Median :  484.8   Median : 0.5784  
##  Mean   : 229.44   Mean   : 1040.95   Mean   :  748.9   Mean   : 1.0430  
##  3rd Qu.: 288.50   3rd Qu.: 1218.88   3rd Qu.:  815.8   3rd Qu.: 1.1648  
##  Max.   :2619.19   Max.   :16100.00   Max.   :29981.4   Max.   :25.9300  
##  NA's   :91        NA's   :50         NA's   :90        NA's   :213      
##       DUR             COLI_FEC             N_NO3                 AS          
##  Min.   :   3.883   Length:2728        Min.   :   0.0310   Min.   :0.000100  
##  1st Qu.: 115.950   Class :character   1st Qu.:   0.7645   1st Qu.:0.001500  
##  Median : 226.551   Mode  :character   Median :   2.0430   Median :0.003461  
##  Mean   : 329.654                      Mean   :  12.7183   Mean   :0.018841  
##  3rd Qu.: 404.787                      3rd Qu.:   4.8335   3rd Qu.:0.018979  
##  Max.   :5855.090                      Max.   :2624.0000   Max.   :1.700000  
##  NA's   :63                            NA's   :105         NA's   :296       
##        CD                CR                HG                PB         
##  Min.   :0.00019   Min.   :0.00100   Min.   :0.00010   Min.   :0.00087  
##  1st Qu.:0.00130   1st Qu.:0.00120   1st Qu.:0.00020   1st Qu.:0.00154  
##  Median :0.00130   Median :0.00140   Median :0.00020   Median :0.00154  
##  Mean   :0.00147   Mean   :0.00672   Mean   :0.00028   Mean   :0.00384  
##  3rd Qu.:0.00130   3rd Qu.:0.00353   3rd Qu.:0.00023   3rd Qu.:0.00154  
##  Max.   :0.15398   Max.   :2.14389   Max.   :0.01961   Max.   :0.08747  
##  NA's   :554       NA's   :553       NA's   :498       NA's   :551      
##        MN                FE             SEMAFORO        
##  Min.   :0.00020   Min.   : 0.00410   Length:2728       
##  1st Qu.:0.00044   1st Qu.: 0.02124   Class :character  
##  Median :0.00244   Median : 0.05545   Mode  :character  
##  Mean   :0.06696   Mean   : 0.18585                     
##  3rd Qu.:0.01510   3rd Qu.: 0.13461                     
##  Max.   :3.65678   Max.   :23.67000                     
##  NA's   :472       NA's   :464
# Ver la estructura (para confirmar si los indicadores son numéricos 
# y los semáforos son factores/texto)
str(datos_limpios2)
## tibble [2,728 × 19] (S3: tbl_df/tbl/data.frame)
##  $ ORGANISMO _DE_CUENCA: chr [1:2728] "RÍO BRAVO" "CUENCAS CENTRALES DEL NORTE" "CUENCAS CENTRALES DEL NORTE" "CUENCAS CENTRALES DEL NORTE" ...
##  $ ESTADO              : chr [1:2728] "COAHUILA DE ZARAGOZA" "DURANGO" "DURANGO" "DURANGO" ...
##  $ ACUIFERO            : chr [1:2728] "PRESA LA AMISTAD" "VILLA JUAREZ" "PRINCIPAL-REGION LAGUNERA" "PRINCIPAL-REGION LAGUNERA" ...
##  $ SUBTIPO             : chr [1:2728] "BROTE" "POZO" "POZO" "POZO" ...
##  $ ALC                 : num [1:2728] 140 149 205 229 168 ...
##  $ CONDUCT             : num [1:2728] 884 329 615 602 379 354 637 663 555 560 ...
##  $ SDT                 : num [1:2728] 562 185 378 322 240 ...
##  $ FLUORUROS           : num [1:2728] 0.82 0.52 0.43 0.5 0.566 ...
##  $ DUR                 : num [1:2728] 506 132 230 222 152 ...
##  $ COLI_FEC            : chr [1:2728] "97" "10" "10" "10" ...
##  $ N_NO3               : num [1:2728] 0.591 0.316 3.182 3.843 0.157 ...
##  $ AS                  : num [1:2728] 0.0015 0.0117 0.0128 0.0125 0.0126 ...
##  $ CD                  : num [1:2728] 0.0013 0.0013 0.0013 0.0013 0.0013 ...
##  $ CR                  : num [1:2728] 0.0012 0.0012 0.0023 0.0021 0.00175 0.0013 0.0012 0.0021 0.0012 0.0012 ...
##  $ HG                  : num [1:2728] 0.000201 0.000116 0.000201 0.000201 0.000201 ...
##  $ PB                  : num [1:2728] 0.00154 0.00154 0.00827 0.00154 0.00154 0.00154 0.00154 0.00154 0.00154 0.00154 ...
##  $ MN                  : num [1:2728] 0.032 0.00021 0.01315 0.00111 0.1047 ...
##  $ FE                  : num [1:2728] 0.011 0.00425 0.19075 0.01655 0.0045 ...
##  $ SEMAFORO            : chr [1:2728] "AMARILLO" "VERDE" "VERDE" "VERDE" ...

Al inspeccionar todos los datos podemos ver que tenemos dos inconvenientes para resolver:

  • El indicador de COLI_FEC es del tipo character
  • Existen valores nulos en varios de los indicadores

Indicador del tipo texto

Revisando los datos nos damos cuenta que el motivo por el cual el indicador COLI_FEC es del tipo texto se debe a que hay datos que utilizan el “.” como delimitador decimal y otros que utilizan “,”.

Suele suceder cuando se toman muestras de muchos lugares distintos y se da a lugar que cada uno utilice su forma de registrar los valores

Vamos a cambiar los delimitadores “,” por “.” para unificar criterios y convertir el indicador a numerico

datos_limpios3 <- datos_limpios2 %>%
  mutate(COLI_FEC = {
    # Convertir a texto para manipularlo
    temp <- as.character(COLI_FEC)
    
    # Cambiar comas por puntos
    temp <- str_replace_all(temp, ",", ".")
    
    # Convertir a número final
    as.numeric(temp)
  })

Manejo de datos nulos

En varios de los indicadores encontramos valores nulos por lo que hay que procesarlos para poder trabajar con los datos.

En este caso se definio por reemplazar los valores nulos por 0, ya que analizando esas muestras se determina que no influye el valor del indicador con el dato nulo en el semaforo final, sino que hay otro indicador que determina fuertemente la calidad del agua.

Igualmente podria estudiarse otros tipos de procesamiento de los datos nulos para comprobar si afectan en el resultado final

datos_finales <- datos_limpios3 %>%
  mutate(across(where(is.numeric), ~replace_na(., 0)))

Analisis de Datos Limpios

Registro de Calidad segun Semaforo

# Contar cuántos registros hay por cada color de semáforo
table(datos_finales$SEMAFORO)
## 
## AMARILLO     ROJO    VERDE 
##      459     1031     1238
# Visualizar la distribución


ggplot(datos_finales, aes(x = SEMAFORO, fill = SEMAFORO)) +
  geom_bar() +
  scale_fill_manual(values = c("VERDE" = "#2ecc71", "AMARILLO" = "#f1c40f", "ROJO" = "#e74c3c")) +
  labs(title = "Distribución de Calidad del Agua", x = "Estado", y = "Cantidad") +
  theme_minimal()

Podemos realizar un resumen estadístico en funcion del grafico anterior:

  • Total de muestras: 2,728 sitios.

    • Verde (Buena calidad): 1,238 sitios (45.4%)

    • Rojo (Incumplimiento): 1,031 sitios (37.8%)

    • Amarillo (En riesgo/Aceptable): 459 sitios (16.8%)

Calidad Comprometida: Menos de la mitad de los sitios monitoreados (45.4%) tienen una calidad de agua óptima (“Verde”).

Alerta Sanitaria: Es preocupante que más de un tercio de las muestras (37.8%) estén en Rojo. Esto indica una presencia crítica de contaminantes (como Fluoruros, Arsénico o Coliformes) que superan los límites permisibles para el consumo humano o uso agrícola.

Tendencia Negativa: Si sumamos las categorías Amarillo y Rojo, tenemos que el 54.6% de los puntos de extracción presentan algún grado de contaminación o sospecha. Esto sugiere que la mayoría del agua subterránea analizada requiere algún tipo de tratamiento antes de ser considerada segura.

Analisis de indicadores por estado

Ahora buscamos encontrar patrones en los indicadores segun el lugar donde se obtuvieron las muestras. Esto nos va ayudar a sacar conclusiones sobre la calidad del agua en los diferentes estados de Mexico, donde vamos a poder ver:

  • Columnas Verticales: Si una columna está muy roja en casi todos los estados, significa que ese contaminante es un problema sistémico en todo México.

  • Filas Horizontales: Si un Estado tiene varios cuadros rojos o naranjas, es una zona crítica donde convergen varios problemas de calidad de agua.

Para mejorar la visualizacion del grafico vamos a ordenar los estados desde el Norte al Sur asi podemos sacar distintos tipos de concluciones rapidas:

  • Correlación Geológica: Podrás observar si el bloque superior (Norte) muestra colores más intensos en Fluoruros o SDT debido a la aridez y la composición del suelo.

  • Patrones de Saneamiento: Podrás ver si el bloque inferior (Sur/Sureste) tiene una mayor intensidad en E. Coli debido a las características de las cuencas hidrológicas o gestión de residuos.

  • Lectura Intuitiva: El lector del gráfico ya no tiene que buscar un estado al azar; su cerebro procesa la información como si estuviera viendo un mapa de México “estirado”.

# Lista de estados ordenados por latitud aproximada (Norte -> Sur)
orden_geografico <- c(
  "BAJA CALIFORNIA", "SONORA", "CHIHUAHUA", "COAHUILA DE ZARAGOZA", "NUEVO LEÓN", "TAMAULIPAS",
  "BAJA CALIFORNIA SUR", "SINALOA", "DURANGO", "ZACATECAS", "SAN LUIS POTOSI", "SAN LUIS POTOSÍ",
  "NAYARIT", "AGUASCALIENTES", "JALISCO", "GUANAJUATO", "QUERÉTARO ARTEAGA", "HIDALGO",
  "COLIMA", "MICHOACÁN DE OCAMPO", "ESTADO DE MÉXICO", "CIUDAD DE MÉXICO", "TLAXCALA",
  "MORELOS", "MORELOS ", "PUEBLA", "VERACRUZ DE IGNACIO DE LA LLAVE", "GUERRERO", "OAXACA", "TABASCO", "CHIAPAS",
  "CAMPECHE", "YUCATÁN", "QUINTANA ROO"
)

Como los indicadores tienen unidades muy distintas (el Flúor se mide en unidades y la Conductividad en miles), si los graficamos juntos sin ajustar, la Conductividad “borraría” a los demás. Por eso, vamos a normalizar los datos (escalarlos de 0 a 1) para que el color represente qué tan grave es el nivel en comparación con el resto del país.

# Preparar datos y aplicar el factor ordenado
df_heatmap <- datos_finales %>%
  group_by(ESTADO) %>%
  summarise(
    Fluoruros = mean(FLUORUROS, na.rm = TRUE),
    Conductividad = mean(CONDUCT, na.rm = TRUE),
    SDT = mean(SDT, na.rm = TRUE), # Sólidos Disueltos Totales
    E_Coli = mean(COLI_FEC, na.rm = TRUE),
    Alcalino = mean(ALC, na.rm = TRUE),
    Dureza = mean(DUR, na.rm = TRUE),
    Nitrogenos = mean(N_NO3, na.rm = TRUE),
    Arsenico = mean(AS, na.rm = TRUE),
    Cadmio = mean(CD, na.rm = TRUE),
    Cromo = mean(CR, na.rm = TRUE),
    Mercurio = mean(HG, na.rm = TRUE),
    Plomo = mean(PB, na.rm = TRUE),
    Manganeso = mean(MN, na.rm = TRUE),
    Hierro = mean(FE, na.rm = TRUE),
  ) %>%
  mutate(across(-ESTADO, ~ rescale(.x))) %>%
  pivot_longer(cols = -ESTADO, names_to = "Indicador", values_to = "Nivel_Normalizado") %>%
  # Convertimos a factor con el orden invertido para que el Norte esté arriba
  mutate(ESTADO = factor(ESTADO, levels = rev(orden_geografico)))

# Graficar
p_heatmap <- ggplot(df_heatmap, aes(x = Indicador, y = ESTADO, fill = Nivel_Normalizado)) +
  geom_tile(color = "white") +
  scale_fill_gradient(low = "#fff5f0", high = "#cb181d") +
  labs(
    title = "Gradiente Geográfico de Contaminación (Norte a Sur)",
    x = "Indicador",
    y = "Estado (Sur <- Norte)",
    fill = "Nivel"
  ) +
  theme_minimal() +
  theme(
    axis.text.x = element_text(angle = 45, hjust = 1),
    panel.grid = element_blank()
  )

ggplotly(p_heatmap)

Identificación de Patrones:

  • En general toda el agua de Mexico es bastante Alcalina (la columna de este indicador tiene colores mas intensos)
  • Los estados al Norte tienen mayor presencia de minerales lo que aumenta los indicadores relacionados (Conductividad, Dureza, Solidos Disueltos) por lo que vemos influencia geologica en la calidad del agua
  • En los estados del centro podemos ver indicadores preocupantes como ECOLI, Mercurio y Plomo. Esto indica que puede ser un estado con produccion agricola y poco tratamiento de sus aguas
  • En la capital podemos ver los niveles mas bajos de los indicadores reflejando que es donde se realizan tratamientos exhaustivos

Analisis de indicadores y sus relaciones

Para poder ver las estadisticas descriptivas de los indicadores vamos a realizar graficos de cajas. Al ser varios indicadores (14) y cada uno tiene tres distribuciones (VERDE, AMARILLO, ROJO), creamos una funcion para graficar las cajas y luego solo utilizamos la funcion para el indicador que necesitamos

# Definimos la función maestra
crear_boxplot_metales <- function(columna, nombre_elemento) {
  
  ggplot(datos_finales, aes(x = SEMAFORO, y = .data[[columna]], fill = SEMAFORO)) +
    # Añadimos los puntos (jitter) para ver la densidad real de los datos
    geom_boxplot(outlier.color = "red", alpha = 0.7, outlier.shape = 0.5) + 
    
    # TRANSFORMACIÓN LOGARÍTMICA
    # trans = "log10" comprime los outliers visualmente
    scale_y_log10(labels = scales::comma) + 
    
    # Paleta de colores consistente
    scale_fill_manual(values = c("ROJO" = "#e74c3c", "AMARILLO" = "#f1c40f", "VERDE" = "#2ecc71")) +
    
    # Etiquetas dinámicas
    labs(
      title = paste("Distribución de", nombre_elemento, "(Escala Logarítmica)"),
      subtitle = "La escala logarítmica permite comparar órdenes de magnitud eliminando el sesgo de outliers",
      x = "Semáforo de Calidad",
      y = "Concentración (mg/L) - Eje Log"
    ) +
    theme_minimal() +
    theme(
      legend.position = "none",
      panel.grid.minor = element_blank() # Limpiamos el fondo para mejor visibilidad
    )
}

Por la misma razon que son varios graficos a realizar, lo hacemos en tandas de a tres para que sean legibles y no ocupen demasiada memoria en el codigo de ejecutar todos.

Vamos a utilizar la escala logaritmica para los graficos ya que los indicadores tienen muchos outliers que hacen incomprensible los grafico en condiciones normales. Al usar la escala logaritmica nos permite:

  • Compresión inteligente: El logaritmo “encoge” los valores gigantes (outliers) pero “estira” los valores pequeños. Esto permite que en el mismo gráfico un valor de 0.001 y uno de 1000 sin que el pequeño desaparezca.

  • Interpretación: No se pierde la relación de magnitud. Si un punto está una unidad por encima de otro en escala logarítmica (base 10), significa que es 10 veces más concentrado.

grafico_fluoruros <- crear_boxplot_metales("FLUORUROS","Fluoruros")
grafico_ALC <- crear_boxplot_metales("ALC","Alcalinidad")
grafico_CONDUCT <- crear_boxplot_metales("CONDUCT","Conductividad")




fig1 <- subplot(ggplotly(grafico_ALC), ggplotly(grafico_fluoruros), ggplotly(grafico_CONDUCT), nrows = 1, titleX = TRUE, titleY = TRUE, margin = 0.07) %>%
  layout(title = "Comparativa de Contaminantes por Semáforo de Calidad", showlegend = FALSE)

fig1
grafico_SDT <- crear_boxplot_metales("SDT","Solidos Disueltos")

grafico_DUR <- crear_boxplot_metales("DUR","Dureza")

grafico_COLI <- crear_boxplot_metales("COLI_FEC","Coliformes")

fig2 <- subplot( ggplotly(grafico_SDT), ggplotly(grafico_DUR), ggplotly(grafico_COLI), nrows = 1, titleX = TRUE, titleY = TRUE, margin = 0.07) %>%
  layout(title = "Comparativa de Contaminantes por Semáforo de Calidad", showlegend = FALSE)

fig2
grafico_NNo3 <- crear_boxplot_metales("N_NO3","Nitratos")

grafico_AS <- crear_boxplot_metales("AS","Arsenico")

grafico_CD <- crear_boxplot_metales("CD","Cadmio")

fig3 <- subplot(ggplotly(grafico_NNo3), ggplotly(grafico_AS), ggplotly(grafico_CD), nrows = 1, titleX = TRUE, titleY = TRUE, margin = 0.07) %>%
  layout(title = "Comparativa de Contaminantes por Semáforo de Calidad", showlegend = FALSE)

fig3
grafico_CR <- crear_boxplot_metales("CR","Cromo")

grafico_HG <- crear_boxplot_metales("HG","Mercurio")

grafico_PB <- crear_boxplot_metales("PB","Plomo")



fig4 <- subplot(ggplotly(grafico_CR), ggplotly(grafico_HG), ggplotly(grafico_PB), nrows = 1, titleX = TRUE, titleY = TRUE, margin = 0.07) %>%
  layout(title = "Comparativa de Contaminantes por Semáforo de Calidad", showlegend = FALSE)

fig4
grafico_MN <- crear_boxplot_metales("MN","Manganeso")

grafico_FE <- crear_boxplot_metales("FE","Hierro")

fig5 <- subplot(ggplotly(grafico_MN), ggplotly(grafico_FE), nrows = 1, titleX = TRUE, titleY = TRUE, margin = 0.07) %>%
  layout(title = "Comparativa de Contaminantes por Semáforo de Calidad", showlegend = FALSE)

fig5

Por los graficos podemos ver la cantidad de valores extremos que tienen todos los indicadores. Esto se va a tener en cuenta en los siguientes analisis para ejecutar metodos que no sigan una distribucion normal

Matriz de correlacion

La matriz de correlacion convencional es la de Pearson pero para este caso vamos a utilizar la de Spearman. Esto se debe a:

1. Correlación de Pearson (El estándar lineal)

  • Mide la fuerza y la dirección de una relación lineal entre dos variables.

  • Qué busca: ¿Si dibujo una línea recta, qué tan cerca están los puntos de ella?

  • Supuestos: * Los datos deben seguir una distribución normal (forma de campana).

  • La relación debe ser una línea recta.

  • Sensibilidad: Es muy sensible a los valores atípicos (outliers). Un solo pozo con una contaminación extrema puede “arruinar” el coeficiente y hacer creer que hay una relación donde no la hay (o viceversa).

2. Correlación de Spearman (El estándar de rangos)

  • Mide la fuerza y la dirección de una relación monótona.

  • Qué busca: ¿Si una variable aumenta, la otra también aumenta (aunque no sea en línea recta)?

  • Cómo funciona: En lugar de usar los valores reales (0.005, 1500, 22), R convierte todo a posiciones (ranks): 1º, 2º, 3º…

  • Supuestos: No asume una distribución normal (no paramétrica).

  • Robustez: Es muy resistente a los valores atípicos. Como usa rangos, a Spearman no le importa si un valor es 1,000 o 1,000,000; simplemente es el valor “más alto” (ej. el puesto 100).

# 2. Seleccionar solo las columnas numéricas de interés
datos_num <- datos_finales %>%
  select(FLUORUROS, CONDUCT, SDT, COLI_FEC, AS, N_NO3, ALC, DUR, CD, CR, HG, PB, MN, FE) # Añade los que tengas

# 3. Calcular la matriz de Spearman (robusta a valores atípicos)
matriz_cor <- cor(datos_num, use = "complete.obs", method = "spearman")

# 4. Graficar la matriz
ggcorrplot(matriz_cor, 
           hc.order = TRUE,         # Agrupa variables similares
           type = "lower",          # Solo muestra la mitad inferior (evita repetición)
           lab = TRUE,              # Muestra el número de la correlación
           lab_size = 3,
           title = "Interacción entre Indicadores de Calidad del Agua",
           colors = c("#6D9EC1", "white", "#E46726"), # Azul (negativa), Blanco (nula), Naranja (positiva)
           outline.color = "white")

Indicadores mas influyentes

Factor de Impacto Relativo

Compararemos el valor promedio de cada contaminante en los sitios Rojos contra el promedio en los sitios Verdes. Si un indicador tiene un promedio 10 veces mayor en las zonas rojas que en las verdes, ese es un “culpable” principal de la degradación del agua.

# Calcular el impacto: Relación de promedios Rojo vs Verde
analisis_impacto <- datos_finales %>%
  group_by(SEMAFORO) %>%
  summarise(across(where(is.numeric), mean, na.rm = TRUE)) %>%
  # Pasamos a formato largo para calcular la proporción
  pivot_longer(cols = -SEMAFORO, names_to = "Indicador", values_to = "Valor") %>%
  pivot_wider(names_from = SEMAFORO, values_from = Valor) %>%
  # Calculamos cuántas veces es más alta la concentración en ROJO que en VERDE
  mutate(Factor_Impacto = ROJO / VERDE) %>%
  # Filtramos indicadores con impacto significativo y ordenamos
  filter(!is.na(Factor_Impacto)) %>%
  arrange(desc(Factor_Impacto))

# Graficar el Impacto
grafico_impacto <- ggplot(analisis_impacto, aes(x = reorder(Indicador, Factor_Impacto), 
                                                y = Factor_Impacto, 
                                                fill = Factor_Impacto)) +
  geom_col() +
  coord_flip() + # Barra horizontal para leer mejor los nombres
  scale_fill_gradient(low = "#f39c12", high = "#c0392b") +
  labs(
    title = "Culpables de la Calidad Roja",
    subtitle = "Factor de incremento (Promedio Rojo / Promedio Verde)",
    x = "Indicador de Contaminación",
    y = "Veces más concentrado en zonas Rojas"
  ) +
  theme_minimal() +
  theme(legend.position = "none")

# Hacerlo interactivo para ver el multiplicador exacto
ggplotly(grafico_impacto)

Resultados:

  • El Multiplicador: Si el indicador COLI_FECAL tiene un factor de 14, significa que en promedio, los pozos en “Rojo” tienen 14 veces más bacterias que los pozos en “Verde”. Esto lo convierte en un factor crítico de salud pública.

  • Químicos vs. Bacteriológicos: Los Fluoruros y el Arsénico tienen los factores más altos por lo que el problema es de origen geológico (el agua está disolviendo minerales del subsuelo).

  • Los Nitratos (NNO3) y Coliformes estando arriba indican que el problema es antropogénico (uso de fertilizantes y falta de tratamiento de aguas residuales).

  • Indicadores de Salinidad: La Conductividad y los SDT (Sólidos Disueltos Totales) suelen tener factores altos cuando hay sobreexplotación de acuíferos, lo que causa intrusión salina.

Importancia de las variables (prueba de Krsukal-Wallis)

Utilizaremos una prueba estadística llamada Kruskal-Wallis, que nos va ayudar porque:

  • No asume que los datos siguen una campana de Gauss.

  • Compara si las concentraciones son realmente distintas entre el Verde, Amarillo y Rojo.

La prueba de Kruskal-Wallis es una prueba no paramétrica que compara las medianas de más de dos grupos independientes, siendo una alternativa al ANOVA cuando no se cumplen los supuestos de normalidad, y evalúa si los datos provienen de la misma población o distribución. Funciona asignando rangos a todos los datos combinados y comparando las sumas de rangos entre los grupos para determinar si hay una diferencia significativa.

# Aseguramos que el Semáforo sea un Factor
datos_finales$SEMAFORO <- as.factor(datos_finales$SEMAFORO)

# Calculamos la importancia para cada indicador
importancia_vars <- datos_finales %>%
  select(where(is.numeric), SEMAFORO) %>%
  pivot_longer(cols = -SEMAFORO, names_to = "Indicador", values_to = "Valor") %>%
  group_by(Indicador) %>%
  summarise(
    # El estadístico chi-squared indica qué tan diferentes son los grupos
    Estadistico_Importancia = kruskal.test(Valor ~ SEMAFORO)$statistic,
    P_Value = kruskal.test(Valor ~ SEMAFORO)$p.value
  ) %>%
  # Ordenamos de mayor a menor importancia
  arrange(desc(Estadistico_Importancia))

# Mostramos el Top 10 de indicadores más influyentes
print(importancia_vars)
## # A tibble: 14 × 3
##    Indicador Estadistico_Importancia   P_Value
##    <chr>                       <dbl>     <dbl>
##  1 FLUORUROS                  727.   1.26e-158
##  2 AS                         562.   1.06e-122
##  3 CONDUCT                    497.   1.07e-108
##  4 SDT                        472.   3.65e-103
##  5 DUR                        369.   8.13e- 81
##  6 FE                         240.   8.55e- 53
##  7 MN                         218.   5.06e- 48
##  8 ALC                        105.   1.83e- 23
##  9 N_NO3                       62.7  2.40e- 14
## 10 PB                          54.1  1.76e- 12
## 11 CD                          41.9  7.84e- 10
## 12 CR                          34.1  3.89e-  8
## 13 HG                          28.3  7.28e-  7
## 14 COLI_FEC                     9.80 7.46e-  3
# Graficar el Ranking
ggplot(importancia_vars, aes(x = reorder(Indicador, Estadistico_Importancia), y = Estadistico_Importancia)) +
  geom_bar(stat = "identity", fill = "#3498db") +
  coord_flip() +
  theme_minimal() +
  labs(title = "Ranking de Indicadores: ¿Quién define el color del agua?",
       subtitle = "Basado en la prueba de Kruskal-Wallis (Mayor valor = Más influyente)",
       x = "Indicador Químico",
       y = "Puntaje de Importancia")

Resultados:

  • Estadistico_Importancia (Chi-squared): Cuanto más alto sea este número, más “poder” tiene ese indicador para cambiar el semáforo. Si un indicador tiene un puntaje de 500 y otro de 10, el de 500 es el que realmente está dictando las reglas.

  • p_Value: Si es menor a 0.05, significa que las diferencias que viste en los boxplots son estadísticamente reales y no fruto de la casualidad.

Vemos que hay diferencias entre una prueba y la otra al determinar que indicador es mas influyente para determinar la calidad del agua. Esto se da porque los valores atípicos (outliers) y la distribución de los datos afectan nuestros resultados.

Diferencias

1. Método del Impacto (Promedio ROJO / Promedio VERDE)

Este método mide la magnitud relativa. Básicamente estás preguntando: “¿Cuántas veces es más grande la concentración en los sitios contaminados comparado con los limpios?

  • Ecoli y Nitratos salen primero: Estos indicadores suelen tener comportamientos “explosivos”. En un sitio con fuga de drenaje, el Ecoli puede pasar de 10 a 100,000 unidades de un momento a otro. Aunque solo ocurra en pocos sitios, esos valores gigantescos “inflan” el promedio del grupo ROJO.

  • Sensibilidad: Es extremadamente sensible a outliers. Un solo pozo terriblemente contaminado puede hacer que el promedio suba tanto que parezca que el indicador es el más influyente de todos, aunque en el 90% de los otros pozos rojos el valor sea bajo.

  • Limitación: No toma en cuenta la variabilidad ni qué tan común es el problema en toda la muestra.

2. Prueba de Kruskal-Wallis

Esta es una prueba no paramétrica que no usa los valores exactos, sino sus rangos (posiciones). “¿Es la distribución de los valores en el grupo ROJO consistentemente mayor que en el grupo VERDE?

  • Fluoruros y Arsénico salen primero: En estos datos, el Arsénico y los Fluoruros son contaminantes de origen geológico. No son “explosivos” como una bacteria, sino que están presentes de forma más constante y generalizada en los sitios que se clasifican como ROJO.

  • Robustez: Al trabajar con rangos, a la prueba de Kruskal-Wallis no le importa si un valor de Ecoli es 1,000 o 1,000,000; para la prueba, ambos son simplemente “valores altos”. Esto elimina el sesgo de los outliers.

Si el Ecoli tiene un impacto alto en promedio pero un Kruskal-Wallis bajo, significa que pocos sitios están MUY contaminados, pero en la mayoría de los sitios rojos, el Ecoli no es el problema principal. En cambio, si el Fluoruro tiene un KW alto, significa que casi todos los sitios rojos tienen más fluoruro que los verdes.

Modelos de Clasificacion

Para crear los modelos de clasificacion vamos a elegir solamente los valores de los indicadores y es muy importante que la columna SEMAFORO este como factor (que R entienda que son las etiquetas a clasificar)

# Seleccionamos las variables numéricas (Features) y nuestro objetivo (Target)

features_list <- c("ALC", "CONDUCT", "SDT", "FLUORUROS", "DUR", "COLI_FEC", 
                  "N_NO3", "AS", "CD", "CR", "HG", "PB", "MN", "FE")

data_modelo <- datos_finales %>%
  mutate(SEMAFORO = as.factor(SEMAFORO)) %>%
  select(all_of(features_list), SEMAFORO) %>%
  drop_na()

KNN

El primer modelo que vamos a crear esKNN (K-Nearest Neighbors). Es un algoritmo basado en distancia. Para clasificar un nuevo pozo de agua, busca a los \(K\) pozos más parecidos (vecinos) en sus niveles químicos y le asigna el color que tenga la mayoría

# Escalamiento de datos (Normalización)
# Es vital para KNN porque se basa en distancias Euclidianas
preProcValues <- preProcess(data_modelo[, features_list], method = c("center", "scale"))
datos_escalados <- predict(preProcValues, data_modelo)

# Dividir en entrenamiento (80%) y prueba (20%)
set.seed(123) # Para que los resultados sean reproducibles
trainIndex <- createDataPartition(datos_escalados$SEMAFORO, p = .8, list = FALSE)
train_data <- datos_escalados[trainIndex, ]
test_data  <- datos_escalados[-trainIndex, ]

# Ejecutar el modelo KNN
# Probaremos con k = 5 (los 5 vecinos más cercanos)
predicciones <- knn(train = train_data[, features_list], 
                    test = test_data[, features_list], 
                    cl = train_data$SEMAFORO, 
                    k = 5)

# Evaluación del modelo
matriz_confusion <- confusionMatrix(predicciones, as.factor(test_data$SEMAFORO))
print(matriz_confusion)
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction AMARILLO ROJO VERDE
##   AMARILLO       71   24     2
##   ROJO            1  156     4
##   VERDE          19   26   241
## 
## Overall Statistics
##                                           
##                Accuracy : 0.8603          
##                  95% CI : (0.8283, 0.8883)
##     No Information Rate : 0.454           
##     P-Value [Acc > NIR] : < 2.2e-16       
##                                           
##                   Kappa : 0.7744          
##                                           
##  Mcnemar's Test P-Value : 4.761e-11       
## 
## Statistics by Class:
## 
##                      Class: AMARILLO Class: ROJO Class: VERDE
## Sensitivity                   0.7802      0.7573       0.9757
## Specificity                   0.9426      0.9852       0.8485
## Pos Pred Value                0.7320      0.9689       0.8427
## Neg Pred Value                0.9553      0.8695       0.9767
## Prevalence                    0.1673      0.3787       0.4540
## Detection Rate                0.1305      0.2868       0.4430
## Detection Prevalence          0.1783      0.2960       0.5257
## Balanced Accuracy             0.8614      0.8712       0.9121

Resultados:

  • Sensibilidad: el agua se agrupa por similitud, pero es muy delicado. Si no escalamos los datos (normalización), el modelo falla porque una variable grande (como la Conductividad) opacaba a las pequeñas (como el Arsénico).

  • Precisión (~88%): Es buena, pero deja un margen de error importante en las fronteras donde un semáforo cambia de color.

Precision vs K

Variando el hiperparametro K podemos obtener diferentes resultados de precision. Para encontrar el que optimiza nuestro modelo podemos variarlo e ir graficandolo

# Definir rango de K a probar
k_values <- 1:30
accuracy_results <- numeric(length(k_values))

# Bucle para probar cada K
for(i in seq_along(k_values)) {
  set.seed(123)
  pred <- knn(train = train_data[, features_list], 
              test = test_data[, features_list], 
              cl = train_data$SEMAFORO, 
              k = k_values[i])
  
  # Calculamos la precisión (Accuracy)
  accuracy_results[i] <- sum(pred == test_data$SEMAFORO) / nrow(test_data)
}

# Crear DataFrame para graficar
df_plot <- data.frame(K = k_values, Accuracy = accuracy_results)

# Graficar con ggplot2
ggplot(df_plot, aes(x = K, y = Accuracy)) +
  geom_line(color = "steelblue", size = 1) +
  geom_point(color = "darkblue", size = 2) +
  theme_minimal() +
  labs(title = "Optimización de KNN: Accuracy vs Valor de K",
       subtitle = "Determinando el número ideal de vecinos para Calidad del Agua",
       x = "Número de Vecinos (K)",
       y = "Precisión (Accuracy)") +
  scale_x_continuous(breaks = k_values)

## Random Forest

El segundo modelo a crear es Random Forest. Consiste en cientos de “Árboles de Decisión” que trabajan en equipo y cada árbol recibe una parte de los datos y una parte de las variables para dar su opinión. Al final, el bosque toma la decisión por mayoría de votos.

# Preparación de datos (usamos la misma limpieza que antes)
datos_rf <- data_modelo %>%
  mutate(SEMAFORO = as.factor(SEMAFORO)) %>%
  select(all_of(features_list), SEMAFORO) %>%
  drop_na()

# Partición de datos
set.seed(123)
trainIndex <- createDataPartition(datos_rf$SEMAFORO, p = 0.8, list = FALSE)
train_data <- datos_rf[trainIndex, ]
test_data  <- datos_rf[-trainIndex, ]

# Entrenar el Bosque Aleatorio
# ntree = 100 (creamos 100 árboles de decisión)
modelo_rf <- randomForest(SEMAFORO ~ ., 
                          data = train_data, 
                          ntree = 100, 
                          importance = TRUE)

# Predicciones
pred_rf <- predict(modelo_rf, test_data)

# Evaluación
print(confusionMatrix(pred_rf, test_data$SEMAFORO))
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction AMARILLO ROJO VERDE
##   AMARILLO       91    3     2
##   ROJO            0  201     0
##   VERDE           0    2   245
## 
## Overall Statistics
##                                           
##                Accuracy : 0.9871          
##                  95% CI : (0.9737, 0.9948)
##     No Information Rate : 0.454           
##     P-Value [Acc > NIR] : <2e-16          
##                                           
##                   Kappa : 0.9794          
##                                           
##  Mcnemar's Test P-Value : 0.0719          
## 
## Statistics by Class:
## 
##                      Class: AMARILLO Class: ROJO Class: VERDE
## Sensitivity                   1.0000      0.9757       0.9919
## Specificity                   0.9890      1.0000       0.9933
## Pos Pred Value                0.9479      1.0000       0.9919
## Neg Pred Value                1.0000      0.9854       0.9933
## Prevalence                    0.1673      0.3787       0.4540
## Detection Rate                0.1673      0.3695       0.4504
## Detection Prevalence          0.1765      0.3695       0.4540
## Balanced Accuracy             0.9945      0.9879       0.9926
# Ver la importancia de las variables
importance_df <- as.data.frame(importance(modelo_rf))
importance_df$Variable <- rownames(importance_df)

ggplot(importance_df, aes(x = reorder(Variable, MeanDecreaseGini), y = MeanDecreaseGini)) +
  geom_bar(stat = "identity", fill = "darkgreen") +
  coord_flip() +
  theme_minimal() +
  labs(title = "Variables más influyentes en el Semáforo del Agua",
       x = "Parámetro Químico",
       y = "Importancia (Gini Index)")

Resultados:

  • Poder Predictivo (~99%): La calidad del agua se define por umbrales (ej. “Si Arsénico > 0.02, entonces Rojo”). Los árboles son perfectos detectando estos cortes exactos.

  • Revelación (Importancia): Nos permitió ver que el Arsénico y el Flúor son los mas influyentes del dataset; ellos dictan el color del semáforo más que cualquier otra variable.

Conclusiones

Para este proyecto de Calidad del Agua, el Random Forest es el mejor aliado. Al ser un problema basado en normas oficiales (límites máximos permitidos), la estructura de “árbol” de este modelo imita perfectamente la lógica humana y legal de clasificación.

KNN es excelente para entender la “geografía” o similitud química, pero para dar un resultado certero y confiable no supera a Random Forest